logo

This script downloads and maps recent eBird data from a Christmas Bird Count circle. It is useful for helping participants to “scout” birds that may occur in their circle/section. This example is run for the West Hennepin 2023 count circle, but it could easily be adapted to other counts.

Setup

Load packages

# Load libraries
library(rebird)
library(dplyr)
library(tidyr)
library(here)
library(ggplot2)
library(sf)
library(ggmap)
library(mapview)
library(DT)

Set parameters

# Name for run (create output directory)
tag <- "West_Hennepin_20231224"

# Define the center coordinates and circle radius
center_lat <- 45.09468
center_long <- -93.63574
center <- c(center_long, center_lat)
radius_km <- 12.07

# How many days back to pull data?
time_days <- 30
# today <- Sys.Date()
today <- "2023-12-24"

Get recent eBird data

# Create a circular buffer using the center and radius
circle <- st_buffer(st_sfc(st_point(center), crs = 4326), dist = radius_km * 1000)
bbox <- st_bbox(circle)

recent_sigtings_file <- here("Projects", "CBC_Scouting", tag, "CBC_sightings_recent.RDS")

if (file.exists(recent_sigtings_file)) {
  CBC_sightings_recent <- readRDS(recent_sigtings_file)
} else {
  # Get recent species
  sp_list <- ebirdgeo(lat = center_lat, lng = center_long, dist = radius_km, back = time_days) %>%
    pull(speciesCode)

  CBC_sightings_recent_raw <- mapply(ebirdgeo, species = sp_list, lat = center_lat, lng = center_long, dist = radius_km, back = time_days)

  # Combine in one dataframe
  CBC_sightings_recent <- bind_rows(
    lapply(CBC_sightings_recent_raw, function(df) {
      if ("exoticCategory" %in% colnames(df)) {
        # If "exoticCategory" column is present, keep it, else create it with NA values
        df <- df %>%
          mutate(exoticCategory = ifelse(!("exoticCategory" %in% colnames(df)), NA, as.character(exoticCategory)))
      } else {
        # If "exoticCategory" column is not present, create it with NA values
        df$exoticCategory <- NA_character_
      }
      return(df)
    })
  )

  saveRDS(CBC_sightings_recent, recent_sigtings_file)
}

Map observations

# Map will bin observations into which week they were observed.
# First project and format date field
CBC_sightings_recent <- CBC_sightings_recent %>%
  st_as_sf(coords = c("lng", "lat"), crs = 4326) %>%
  mutate(Date = as.Date(obsDt))

# Bin dates (approximate number of weeks)
num_bins <- round(time_days / 7, 0)

# Calculate bin edges, extending the range slightly beyond the minimum and maximum dates
min_date <- min(CBC_sightings_recent$Date) - 1
max_date <- max(CBC_sightings_recent$Date)
bin_edges <- seq(min_date, max_date, length.out = num_bins + 1)

# Bin the Date variable
CBC_sightings_recent$DateBin <- cut(CBC_sightings_recent$Date, breaks = bin_edges, labels = FALSE, right = TRUE)

# Get a map tile from a map service (e.g., Stadia Alidade). This ill throw an error if API key isn't in .Renviron.
map_tile <- get_stadiamap(bbox = c(bbox[[1]], bbox[[2]], bbox[[3]], bbox[[4]]), zoom = 11, maptype = "alidade_smooth")

# Generate map
p <- ggmap(map_tile) +
  geom_sf(data = circle, fill = "transparent", color = "grey20", size = 1, inherit.aes = FALSE) +
  geom_sf(data = st_jitter(CBC_sightings_recent), size = 1.5, alpha = .9, aes(color = as.factor(DateBin)), inherit.aes = FALSE) +
  scale_color_manual(name = "Week starting:", values = rev(viridisLite::inferno(num_bins)), labels = format(bin_edges, "%Y-%m-%d")) +
  ggthemes::theme_map() +
  theme(
    legend.position = "top",
    legend.justification = "right",
    legend.background = element_rect(
      fill = "grey90",
      linewidth = 0.25, linetype = "solid",
      colour = "grey40"
    ),
    legend.key = element_rect(fill = "grey90"),
    plot.title = element_text(size = 8, face = "bold"),
    plot.subtitle = element_text(size = 6, face = "italic"),
    strip.text = element_text(size = 6)
  ) +
  facet_wrap(~comName) +
  labs(
    title = "West Hennepin CBC: recent eBird observations",
    subtitle = paste0("Data from 30 day period ending ", today)
  )

ggsave(plot = p, here("Projects", "CBC_Scouting", tag, paste0(tag, "_recent_ebird_data.pdf")), width = 8.5, height = 11)

p

Generate hotspot summaries

recent_hotspots_sigtings_file <- here("Projects", "CBC_Scouting", tag, "CBC_hotspot_sightings_recent.RDS")
hotspots_file <- here("Projects", "CBC_Scouting", tag, "CBC_hotspots.RDS")

if (file.exists(recent_hotspots_sigtings_file) & file.exists(hotspots_file)) {
  CBC_hotspot_sightings_recent <- readRDS(recent_hotspots_sigtings_file)
  CBC_hotspots <- readRDS(hotspots_file)
} else {
  CBC_hotspots <- rebird::ebirdhotspotlist(lat = center_lat, lng = center_long, dist = radius_km) %>%
    st_as_sf(coords = c("lng", "lat"), crs = 4326)
  saveRDS(CBC_hotspots, hotspots_file)

  mapview(CBC_hotspots, label = "locName")

  CBC_hotspot_sightings_recent <- sapply(CBC_hotspots$locId, ebirdregion, simple = FALSE, back = time_days)

  # Combine in one data frame
  CBC_hotspot_sightings_recent <- bind_rows(
    lapply(CBC_hotspot_sightings_recent, function(df) {
      if ("exoticCategory" %in% colnames(df)) {
        # If "exoticCategory" column is present, keep it, else create it with NA values
        df <- df %>%
          mutate(exoticCategory = ifelse(!("exoticCategory" %in% colnames(df)), NA, as.character(exoticCategory)))
      } else {
        # If "exoticCategory" column is not present, create it with NA values
        df$exoticCategory <- NA_character_
      }
      return(df)
    })
  )

  CBC_hotspot_sightings_recent <- CBC_hotspot_sightings_recent %>%
    replace_na(list(howMany = 1))

  saveRDS(CBC_hotspot_sightings_recent, recent_hotspots_sigtings_file)
}

# Generate table with recent species observed at each hotspot (with date last documented)
species_by_hotspot <- CBC_hotspot_sightings_recent %>%
  mutate(obsDtmd = format(as.Date(obsDt), "%m-%d")) %>%
  pivot_wider(id_cols = comName, names_from = locName, values_from = obsDtmd, values_fn = max) %>%
  left_join(rebird:::tax %>% select(comName, taxonOrder)) %>%
  arrange(taxonOrder) %>%
  select(-taxonOrder)
species_by_hotspot[is.na(species_by_hotspot)] <- ""
species_by_hotspot <- species_by_hotspot %>%
  select(sort(tidyselect::peek_vars())) %>%
  select(comName, everything())

species_by_hotspot %>%
  datatable(
    options = list(pageLength = 25),
    caption = "Hotspots each species has been recently documented in"
  )
write.csv(species_by_hotspot, here("Projects", "CBC_Scouting", tag, "species_by_hotspot_last_seen.csv"))

count_hotspots_by_sp <- CBC_hotspot_sightings_recent %>%
  group_by(comName) %>%
  summarize(
    n_hotspots_observed = length(unique(locName)),
    parks_observed = paste(sort(unique(locName)), collapse = ", ")
  ) %>%
  left_join(rebird:::tax %>% select(comName, taxonOrder)) %>%
  arrange(taxonOrder) %>%
  select(-taxonOrder) %>%
  arrange(desc(n_hotspots_observed))

count_hotspots_by_sp %>%
  datatable(
    options = list(pageLength = 25),
    caption = "Number of hotspots each species has recently been observed in"
  )
write.csv(species_by_hotspot, here("Projects", "CBC_Scouting", tag, "count_hotspots_by_sp_raw.csv"))

# rare sp are those seen in fewer than 30% of hotspots
rare_sp <- count_hotspots_by_sp %>%
  filter(n_hotspots_observed <= nrow(CBC_hotspots) * .3) %>%
  pull(comName)

# unique sp are those seen in only 1 hotspot
unique_sp <- count_hotspots_by_sp %>%
  filter(n_hotspots_observed == 1) %>%
  pull(comName)

count_sp_by_hotspot <- CBC_hotspot_sightings_recent %>%
  group_by(locName) %>%
  summarize(
    n_sp = length(unique(comName)),
    n_rare_sp = length(intersect(rare_sp, unique(comName))),
    n_unique_sp = length(intersect(unique_sp, unique(comName))),
    rare_sp = paste(sort(intersect(rare_sp, unique(comName))), collapse = ", "),
    unique_sp = paste(sort(intersect(unique_sp, unique(comName))), collapse = ", ")
  ) %>%
  arrange(desc(n_sp))

write.csv(species_by_hotspot, here("Projects", "CBC_Scouting", tag, "count_sp_by_hotspot_raw.csv"))

count_sp_by_hotspot %>%
  datatable(
    options = list(pageLength = 25),
    caption = "Hotspots with total number of recently species and rare/unique species. Rare species are those seen in <=30% of hotspots."
  )

Session info

sessionInfo()
## R version 4.2.2 (2022-10-31 ucrt)
## Platform: x86_64-w64-mingw32/x64 (64-bit)
## Running under: Windows 10 x64 (build 22621)
## 
## Matrix products: default
## 
## locale:
## [1] LC_COLLATE=English_United States.utf8  LC_CTYPE=English_United States.utf8   
## [3] LC_MONETARY=English_United States.utf8 LC_NUMERIC=C                          
## [5] LC_TIME=English_United States.utf8    
## 
## attached base packages:
## [1] stats     graphics  grDevices utils     datasets  methods   base     
## 
## other attached packages:
##  [1] mapview_2.11.2          ggmap_4.0.0             rebird_1.3.0            DT_0.31                
##  [5] CoordinateCleaner_3.0.1 rgbif_3.7.8             MASS_7.3-58.1           viridis_0.6.4          
##  [9] viridisLite_0.4.2       ggrepel_0.9.4           auk_0.7.0               readxl_1.4.3           
## [13] kableExtra_1.3.4        forcats_1.0.0           stringr_1.5.1           purrr_1.0.2            
## [17] readr_2.1.4             tidyr_1.3.0             tibble_3.2.1            tidyverse_2.0.0        
## [21] lubridate_1.9.3         sf_1.0-14               rgdal_1.6-4             sp_2.1-2               
## [25] RODBC_1.3-23            effects_4.2-2           carData_3.0-5           ggthemes_5.0.0         
## [29] ggplot2_3.4.4           dplyr_1.1.4             knitr_1.45              here_1.0.1             
## 
## loaded via a namespace (and not attached):
##   [1] systemfonts_1.0.5   plyr_1.8.9          lazyeval_0.2.2      splines_4.2.2       crosstalk_1.2.1    
##   [6] leaflet_2.2.1       urltools_1.7.3      digest_0.6.33       htmltools_0.5.7     fansi_1.0.5        
##  [11] magrittr_2.0.3      tzdb_0.4.0          vroom_1.6.4         svglite_2.1.2       timechange_0.2.0   
##  [16] jpeg_0.1-10         colorspace_2.1-0    rvest_1.0.3         mitools_2.4         textshaping_0.3.7  
##  [21] leafem_0.2.3        xfun_0.41           crayon_1.5.2        jsonlite_1.8.7      lme4_1.1-35.1      
##  [26] survival_3.4-0      glue_1.6.2          gtable_0.3.4        webshot_0.5.5       scales_1.3.0       
##  [31] oai_0.4.0           DBI_1.1.3           Rcpp_1.0.11         units_0.8-5         bit_4.0.5          
##  [36] proxy_0.4-27        stats4_4.2.2        survey_4.2-1        htmlwidgets_1.6.4   epuRate_0.1        
##  [41] httr_1.4.7          geosphere_1.5-18    wk_0.9.1            ellipsis_0.3.2      pkgconfig_2.0.3    
##  [46] farver_2.1.1        nnet_7.3-18         sass_0.4.7          utf8_1.2.4          crul_1.4.0         
##  [51] tidyselect_1.2.0    labeling_0.4.3      rlang_1.1.2         munsell_0.5.0       cellranger_1.1.0   
##  [56] tools_4.2.2         cachem_1.0.8        cli_3.6.1           generics_0.1.3      evaluate_0.23      
##  [61] fastmap_1.1.1       ragg_1.2.6          yaml_2.3.7          bit64_4.0.5         s2_1.1.4           
##  [66] satellite_1.0.4     nlme_3.1-160        whisker_0.4.1       xml2_1.3.5          compiler_4.2.2     
##  [71] rstudioapi_0.15.0   curl_5.1.0          png_0.1-8           e1071_1.7-13        bslib_0.6.1        
##  [76] stringi_1.8.2       highr_0.10          lattice_0.20-45     Matrix_1.6-4        classInt_0.4-10    
##  [81] nloptr_2.0.3        vctrs_0.6.5         pillar_1.9.0        lifecycle_1.0.4     triebeard_0.4.1    
##  [86] jquerylib_0.1.4     data.table_1.14.8   bitops_1.0-7        insight_0.19.7      raster_3.6-26      
##  [91] R6_2.5.1            KernSmooth_2.23-20  gridExtra_2.3       codetools_0.2-18    boot_1.3-28        
##  [96] assertthat_0.2.1    rprojroot_2.0.4     withr_2.5.2         httpcode_0.3.0      rnaturalearth_0.3.4
## [101] parallel_4.2.2      hms_1.1.3           terra_1.7-55        grid_4.2.2          class_7.3-20       
## [106] minqa_1.2.6         rmarkdown_2.25      base64enc_0.1-3
 




By Sam Safran